Class 12 Assignment - Mapping Urban Landcover with Spectral Indices and the Surface Temperature band from the Landsat 8 ARD Program

Urban Heat Island Effect source:https://www.epa.gov/heatislands/learn-about-heat-islands

Class 12 Concepts & Themes:

This week’s assignment will encompass the following concepts covered in Class 12 lecture & lab:

  • Application of Spectral Indices to Landsat ARD imagery
  • Utilization of Surface Temperature produce from Landsat ARD
  • Legend Design for Spectral Indices

Class 12 Readings:

-While there are many contextual references in the assignment introduction below, there are no required readings this week, and there will be NO quiz for this week’s materials.

Class 11 Assignment - Preamble:

Assignment Introduction: In the era of climate change, there has been significant research into the earth’s surface temperature, especially in urban areas where socio-economic, landcover, historical urban design all strongly correlate with increasing surface temperatures - often referred to as Urban Heat Islands.

Urban Correlations Between Vegetative Cover & Income

  • Through the historical urban design lens, what is referred to as ‘red lining’ also often correlates with less vegetative cover:

Redlining & Vegetative Cover

Landcover and Surface Temperature in Urban Environments

Landcover Types in Urban Environment Exhibiting Different Surface Temperature Profiles

  • Utilizing the Landsat ARD program, this assignment utilizes first the surface reflectance values across several bands to derive a spectral index - particularly NDVI. To follow, surface temperature values are collected for the same area of interest (AOI), processed to a scaled and offset value, expressed as degrees fahrenheit (°F). The result is a visual comparison of the vegetative cover within the urban AOI relative to the surface temperature result with the expected behavior of correlating low vegetative cover with elevated surface temperature values. Finally, the NDVI index and surface temperature array will be presented using a gradient color map, respective to each theme - vegetative cover and surface temperature in the chosen urban AOI

Assignment Part I - Acquire Landsat 8 ARD Imagery for both surface reflectance and temperature values

  • To Start, this assignment will rely on both the general process and sample imagery found in this week’s lecture demo. Here the city of Springfield, MO is the AOI; the imagery date is the height of summer on June 20th, 2021 at approximately 10:50 am (more on how to determine time of the imagery later in the assignment).

  • For your own assignment, you will pick a US city as your AOI. Good picks are cities that are mid-sized to large. They must first meet the following criteria:

  1. The AOI must fully fit within one ARD grid tile.
  2. Within the ARD grid tile, the city must be covered by just 1 landsat ARD scene.

The Springfield, MO AOI meets both criteria as seen below:

  • Search for a city that falls within one ARD grid tile. Once found, select the grid tile and notate the h and v values in the attribute table for the grid tile. In the case of the grid tile containing Springfield, MO, the h value is 18 and the v value is 11.

ARD Grid overlay to QuickMapServices Google Maps Baselayer

  • Next, Navigate to EarthExplorer, making sure you are logged in as a user (we did this step in Assignment 11 also). Under Search Criteria, Navigate to the Landsat ARD tab for Collection 2 data:

Landsat C2 US ARD

  • Utilize the horizontal and vertical position of the ARD Grid to produce results just for that one grid space. The AOI is Springfield, Missouri which has a Horizontal ID of 18, and Vertical ID of 11:

Horizontal/Vertical Position ID of AOI ARD Tile

  • Next, the following four conditions are met from the returned results:
  1. Cloud-free over Area of Interest - AOI - Springfield, Missouri as example
  2. The tile preview shows full pixel coverage over AOI
  3. The ARD tile date is within ‘the growing season’, i.e. approximately late May through early August. Ideally, the ARD tile date is also in ‘high summer’ so that the surface temperature values will represent a typical hot day at the AOI.
  4. The ARD scene is part of the Landsat 8 Program, and the ARD scene is part of the Collection 2 segment of the program.

This Tile is cloud-free, full coverage over the AOI, and is within the growing season and the high summer - 6/20/21

  • The Download features two products - the Surface Reflectance Values; and secondly the Surface Temperature values. Make sure to download both bundles, untar them with keka or 7zip and place the uncompressed files into the ARD_SR and ARD_ST, respectively:

Surface Reflectance Bundle

Surface Temperature Bundle

  • **Uncompress into the assignment directory, and create a folder outputs for processing outputs upcoming in the assignment:

Assignment Data Structure

Assignment Part II - Verify Landsat 8 ARD Imagery Date, Time, Azimuth and Elevation

  • As we discussed in Lecture 11, the Landsat program instrument is a passive sensor which requires the sun to effectively record electromagnetic values. As Landsat 8 traverses the continental US, it travels in a north>south direction, and it passes the equator heading southward mid-morning. As a result, the Landsat sensor should pass over a city in the continental US in the earlier part of the day. We can check the exact time that a Landsat scene is recorded, and we can also find the sun’s elevation and azimuth. The following graphic shows both elevation and azimuth where 0° elevation is at sunrise and 90° elevation when the sun is directly overhead:

Azimuth & Elevation Diagram

Landsat data acquisition times are expressed in Greenwich Mean Time (GMT), and can be converted to CST, PST and EST. Since Springfield, MO is CST, we subtract 6 hours from the GMT expression in the metadata:

  • CST = -6 GMT
  • PST = -8 GMT
  • EST = -5 GMT

We can also find the weather of a certain day for a certain city:

For Springfield, MO on the data of imagery acquistion (6-20-21) the temperature profile is as follows (hot by 10am!):

Temperature Profile on Date of Image Acquisition

To Find the Date, Time, Elevation and Azimuth, access the .json file located in either the ARD_ST folder:

.json metadata file for ARD tile scene

  • This file can be opened in a text editor. Once open, look for the section titled “SCENE_METADATA”:

Sections Highlighted in Magenta Include Date, Time, Elevation and Azimuth

From this information, we can place the assignment example for Springfield, MO as follows:

  • The image is taken on 6/20/21 at 10:48 AM.
  • The Sun’s Elevation is relatively ‘high in the sky’ at 67 degrees.
  • The Sun’s position is in the south-eastern sky with a Azimuth of 121 degrees.

In effect, if we were standing facing south in Springfield, MO when this Landsat scene was recorded it would have been a hot morning approach 90 degrees; the sun would have been relatively high in the sky approaching a typical ‘noon’ position at 90 degrees; and the sun’s east>west orientation would be southeast in the sky. Indeed, these conditions are pretty good to represent a typical ‘early, hot summer day’ in the midwest of the US in the city of Springfield, MO. This is the condition that we are seeking for our surface temperature analysis upcoming in the assignment. It is also appropriate for the vegetative cover analysis where we would expect the tree coverage to represent a ‘full-leaf’ summer day.

  • Once the metadata has been checked for appropriateness of imagery to the analysis task, you can proceed with confidence that your Landsat imagery is portraying the date,time, ambient temperature and sun position that we seek: a typical warm to hot summer day with typical, full sun exposure.

Assignment Part III- Processing the Imagery for NDVI

  • To Start, Like the Lecture Demo 12, the NDVI spectral index will be created based on the following formula:

NDVI = (NIR Band - Red Band / (NIR Band + Red Band)

  • As translated to the Landsat 8 band numbers:

NDVI = (Band 5 – Band 4) / (Band 5 + Band 4)

Utilize Bands 4 and 5 for NDVI

  • The following formula utilizes the bands for the SR values in Band 4 and Band 5 located in the ARD_SR folder:

("LC08_CU_018011_20210620_20210703_02_SR_B5@1"-"LC08_CU_018011_20210620_20210703_02_SR_B4@1")/("LC08_CU_018011_20210620_20210703_02_SR_B5@1"+"LC08_CU_018011_20210620_20210703_02_SR_B4@1")

  • Review the results with symbolization applied for values within the NDVI Index:

NDVI low > high Symbolization

Save the results as your_city_name_NDVI.tif into the outputs folder. You can also save a QGIS project into the assignment folder before proceeding to the next steps.

Assignment Part IV- Processing the Imagery for Surface Temperature

  • With NDVI complete, turn to calculating the Surface Temperature utilizing the ST layer in the ARD_ST folder. Import this band layer into QGIS:

ST_B10 Layer

Loaded to QGIS Layers Panel

Surface Temperature (ST) – Represents the temperature of the Earth’s surface in Kelvin (K).

  • Also noted in the guide is a scale factor that is utilized for the product:

Kelvin Values for the ST band

Application of Scale and Offset

  • More information on Scaling HERE

  • ST Band 10 in QGIS, Layers Panel:

Note the unscaled Kelvin Values 28804 - 49411

  • Here the digital numbers need to be scaled accordingly, and an offset of 149 is also applied. If the value is to be expressed in Celcius or Fahrenheit, a further conversion needs to take place via the raster calculator. This will be accomplished in two steps:

    • Scale the raster values in raster calculator

    • Apply the following temperature conversion where K is the scaled raster values and F is the Fahrenheit temperature measurement:

    • F = 1.8(K - 273) + 32

  • To Start, open the Raster Calculator and apply the following scale factor and offset value as noted in the product metadata, exporting the product to the outputsdata folder as st_scaled:

    • "LC08_CU_018011_20210620_20210703_02_ST_B10@1" * 0.00341802 + 149

Scaled + Offset Kelvin Values

  • Next, apply the temperature conversion based on the formula for kelvin to fahrenheit units, and output as your_city_name_ST in the typical, default GeoTiff format.

    • 1.8 * ("st_scaled@1" - 273) + 32

  • Review the resulting raster values, Kelvin vs Fahrenheit:

Kelvin vs Fahrenheit

  • Finally, apply an appropriate color ramp with the city focused in the canvas extent to ensure the surface temperature pattern appears correct where high values are located generally over impervious, high built-up areas; lower values in more suburban and rural areas at often towards the edge of the city:

Utlize Value Tool & QuickMapServices Base Imagery to Evaluate Resulting Surface Temperatures Across the Urban Extent

Assignment Part V - Mask Analysis Layers to a Common AOI Extent

  • At this juncture, both the NDVI and Surface Temperature Analysis Rasters have been created; but there extent is simply the extent of the original Landsat tile edge. As the goal of the assignment is to compare and evaluate the vegetative cover within the respective urban AOI relative to surface temperature, ‘cutting’, i.e. masking to the city AOI is valuable. To do so, we will query OpenStreetMap for the ID for the respective city; acquire this ID in Overpass Turbo and utilize it as polygon mask in QGIS. We’ve done this type of query during the Lecture Demo 10.

  • To Start, navigate to OpenStreetMap, and search for the city by name:

Search by City Name and State

  • Next, record the ID; in the case of Springfield, MO, this is a Relation geometry:

ID of Feature in OSM

  • Navigate to Overpass Turbo and enter the query code for the Feature ID (replace your city ID with that of Springfield, MO in the example):
/*
This query looks for a node, way or relation 
with the given id.
*/
[out:json][timeout:25];
// gather results
(
  // query part for: “id:141244”
  relation(141244);
);
// print results
out body;
>;
out skel qt;

Overpass Turbo Query Result for Springfield, MO

  • Download the result as GeoJSON and place into the outputs folder.

  • Import the GeoJSON file into QGIS:

OSM Springfield, MO Overpass Turbo Query Feature Loaded to QGIS

  • At this juncture, its likely that the boundary feature needs to be fixed before utilizing as a mask feature. Run the export feature through the Fix Geometries tool via the Processing Toolbox and produce a temporary layer:

Fix Geometries

  • Next, mask both the NDVI and Surface Temperature rasters and save them to the outputs folder as NDVI_final.tif and ST_final.tif, respectively:

Process the ST layer, then the NDVI layer to the Fix Geometries Mask using the parameters in the tool

  • Before Proceeding to the final cartographic design, check the two results; they should be both clipped to the city boundary located at OpenStreetMap:

Correctly Masked Layers for NDVI and ST Values

Part VIII - Final Cartography:

  • With both Raster Layers Prepared, save the QGIS project before proceeding. The map layout will feature a side-by-side comparison of the NDVI results to the Surface Temperature results. A Gradient Color Ramp will be provided for each raster theme; both will be displayed with the same map scale, so 1 scale bar can apply to both map layout frames.

  • To Start, like Class 5 assignment, create a Gradient Color Ramp for the NDVI Spectral Index and Surface Temperature Gradient using the following short video as guide:

Pseudocolor Color Ramps Applied to Each Analysis Layer

Note: Make sure to apply a pseudocolor color ramp for the legend development. Like assignment 5, below is a short video outlining the process in QGIS LTR 3.1. In later releases - 3.2x - the default symbolization for rasters is indeed a pseudocolor color ramp, and thus no need to develop a secondary gradient color ramp in the layout.

Gradient Color Ramp

  • With the pseudocolor color ramps applied to each layer, proceed to Map Layout and position both layers into a horizontal, landscape orientation:

Page Properties for Landscape Orientation

  • Next lay in the two analysis layers, and apply a continuous gradient legend to each. Make sure to retain the exact same map scale for each map frame, locking each map frame layer so that the next layer does not draw onto the subsequent layer:

Locking the Layer within the Map Frame

  • Use the following assignment draft as a general template for your map layout:

Assignment 12 Draft Layout Example

Note: in the layout example provided, both legends are oriented with low values at top, high values at bottom. There’s a strong argument to flip these legends so they correspond better. Low NDVI values correspond with high surface temperatures; it would be appropriate and logical to flip the surface temperature legend so that its high values are on top and correspond with the low NDVI values at top of the respective NDVI legend.

The final map design should include the following:

  • Landscape or portrait orientation map sheet (your choice)

  • Main map frame featuring the difference pixels atop quickmapservices basemap imagery.

  • An inset map denoting location of the project area.

  • Titling that explains the location, date range and thematic purpose of the map.

  • Simple legend for the difference class (-1).

  • Source the data as follows: Landsat-8 and Landsat-7 analysis imagery courtesy of the U.S. Geological Survey